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1  Introduction 


Resin  transfer  molding  (RTM)  is  being  increasingly  employed  in  the  man¬ 
ufacture  of  large  components  of  fiber-reinforced  composite  materials.  The 
process  involves  injection  of  resin  into  a  mold  cavity  filled  with  fiber  mat. 
Since  good,  close  tolerances  can  be  obtained  by  a  good  mold,  this  process  is 
suited  for  manufacturing  situations  involving  close  tolerances. 

Most  of  the  present  work  in  RTM  simulations  employs  Eulerian  fixed 
meshes  modeling  the  mold  cavity  with  a  finite  element-control  volume  ap¬ 
proach  and  philosophy  to  track  the  flow  field  and  filled  regions  of  the  mold 
[1-5].  In  these  finite  element-control  volume  based  approaches,  a  parameter 
called  fill  factor  is  associated  with  each  node.  Each  node  is  associated  with 
a  control  volume  region,  and  the  nodal  fill  factor  indicates  if  the  associated 
control  volume  region  is  empty  or  full.  A  fill  factor  of  1  corresponds  to  a 
filled  region  saturated  with  the  resin,  while  a  fill  factor  of  0  corresponds  to 
an  unfilled  region.  The  resin  flow  front  exists  in  control  volume  regions  where 
the  fill  factor  is  between  0  and  1.  The  solution  process  involves  solving  for  the 
pressure  field  based  on  both  the  continuity  equation  invoking  Darcy’s  law  and 
the  fill  factor.  The  filled  regions  and  the  fill  factors  are  computed  based  on 
the  actual  flow  rates  and  the  control  volumes  associated  with  each  node.  This 
approach  makes  the  mold-filling  process  to  be  regarded  as  a  quasi-steady  pro¬ 
cess,  even  though  it  is  a  transient  process,  assuming  a  steady-state  condition 
at  each  time  step.  At  each  time  step,  the  pressure  distribution  is  computed 
based  on  both  the  boundary  conditions  and  the  fill  factor  condition.  The 
velocity  field,  hence,  the  flow  rate  computed  from  the  pressure  field,  is  used 
to  update  the  new  front  location,  fill  stage,  and  time  of  simulation.  The  se¬ 
lection  of  the  time  step  increment  for  each  of  the  quasi-steady  state  is  based 
on  the  consideration  that  the  time  increment  used  allows  only  one  control 
volume  region  to  be  completely  filled.  This  restriction  of  the  time  increment 
ensures  the  stability  of  the  quasi-steady  state  approximation  and  forms  the 
so-called  explicit  filling  methodology.  In  principle,  three  approaches  are  evi¬ 
dent  in  the  filling  simulations  based  on  Eulerian  meshes,  which  includes  the 
finite-element-control  volume  method  and  other  related  techniques  such  as 
the  marker  and  cell  (MAC)  method  and  the  volume  of  fluid  (VOF)  method. 

Other  techniques  to  determine  the  filled  regions  have  been  employed  in 
analyzing  the  flow  and  filled  regions  in  metal  casting  flow  simulations  [6]. 
The  MAC  technique  combines  the  computational  Eulerian  mesh  with  a  set 
of  Lagrangian  particles  (the  markers)  representing  the  fluid  configuration, 
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and  moving  through  the  Eulerian  mesh  according  to  the  velocity  vector  field 
and  defining  where  the  fluid  is  located.  In  these  metal  casting  flow  simula¬ 
tions,  the  transient  Navier-Stokes  and  continuity  equations  are  solved  for  the 
velocity  field.  The  VOF  is  another  method  that  is  commonly  employed  in 
mold  flow  simulation  in  metal  castings[7-10].  In  this,  a  factor,  F,  represent¬ 
ing  the  fill  factors  associated  with  the  node  and  the  control  volume  region, 
is  tracked  and  solved  using  a  transport  equation  of  the  form: 

f  +  -VF  =  0.  (1) 

The  velocity  field  is  again  solved  using  the  time-dependent  Navier-Stokes 
equation  and  the  continuity  equation,  and  the  factor  F  is  between  0  and  1. 

In  considering  the  aforementioned  approaches  employed  in  the  filling  of 
mold  cavities  either  in  polymer-filling  or  casting  applications,  it  is  clear  that, 
in  addition  to  the  flow  field,  a  fill  factor  indicating  the  status  of  the  region 
(filled  or  empty)  is  to  be  solved  and  determined.  The  explicit  filling  method¬ 
ology,  described  earlier  and  used  in  finite  element-control  volume  techniques 
for  resin  filling  inside  a  porous  mold  cavity,  involves  only  the  continuity  equa¬ 
tion  given  by  =  0,  based  on  an  incompressible  resin  and  the  velocity  field 
approximated  using  Darcy’s  approximation.  This  necessitates  the  solution  of 
the  pressure  field  using  the  continuity  equation  and  explicit  determination  of 
velocity  fields  and  the  flowrates  based  on  the  velocity  field  in  each  of  the  con¬ 
trol  volume  regions  defined  for  each  node  of  the  finite  element  mesh  defined 
for  the  pressure  field  solution.  This  approach  also  restricts  the  incremental 
time  employed  in  the  time  marching  of  the  resin  flow  front. 

It  is  to  be  noted  that  the  objective  of  the  mold-filling  analysis  is  the 
conservation  of  resin  mass  at  any  instant  of  time  and  determination  of  its 
distribution  inside  a  Eulerian  mesh  mold  cavity.  The  fill  factor  implicitly 
defines  the  amount  of  resin  present  at  any  time  inside  a  mold  cavity  and 
its  distribution.  In  this  report,  an  alternate  form  of  the  mass  balance,  or 
continuity  equation,  involving  a  time  derivative  of  fill  factor,  to  be  determined 
is  proposed.  The  velocity  field  is  then  approximated  based  on  Darcy’s  law  for 
flow  through  a  porous,  permeable  medium.  Finite  element  approximations 
for  both  the  pressure  field  and  the  fill  factor  are  introduced,  thus  providing 
a  pure  finite  element  methodology  not  involving  the  definition  of  the  control 
volumes  for  both  the  pressure  field  and  the  fill  factors.  The  resulting  time- 
dependent  system  of  equations  are  solved  both  for  the  pressure  field  and  the 
fill  factors  defining  the  resin  location  inside  a  Eulerian  mold  cavity  in  an 
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iterative  manner.  The  time-dependent  system  of  equations  can  be  solved 
employing  any  time  step,  and  the  pressure  field  and  the  filled  regions  at 
discrete  time  steps  during  the  analysis  are  determined. 

The  developments  are  presented  based  on  the  following  assumptions:  1) 
the  polymer  resin  is  incompressible;  2)  the  viscosity  of  the  resin  is  taken  to 
be  constant  and  the  filling  is  taken  to  be  at  isothermal  conditions;  3)  the 
flow  is  governed  by  Darcy’s  law  and  the  Reynolds  number  is  small  to  neglect 
the  effect  of  inertial  terms;  and  4)  the  body  forces  and  the  gravity  effects 
are  not  included.  Though  the  developments  are  presented  based  on  these 
assumptions,  the  present  formulations  can  be  easily  extended  to  situations 
involving  nonisothermal  cases  and  where  the  body  forces  and  gravity  effects 
are  involved. 

The  accuracy  issues  and  computational  advantages  of  the  pure  finite- 
element  methodology  proposed  here  are  highlighted.  Several  numerical  ex¬ 
amples  are  presented  to  demonstrate  the  applicability  of  the  method  for 
polymer-filling  simulations.  Though  the  examples  are  presented  in  the  con¬ 
text  of  thin  mold  configurations  in  general  3-D  space,  the  developments  are 
valid  for  general  1-,  2-,  or  3-D  mold-filling  simulations  and  will  be  demon¬ 
strated  in  a  subsequent  report. 


2  Pure  Finite  Element  Methodology 

Before  presenting  the  pure  finite  element  methodology  for  the  solution  of 
the  pressure  field  and  the  fill  factors,  an  alternate  form  of  conservation  mass 
equation  is  considered.  This  equation  is  based  on  resin  mass  conservation 
involving  the  time  derivative  of  the  fill  factor. 

2.1  Conservation  of  Resin  Mass 

As  discussed  earlier,  the  objective  of  the  polymer  resin  mold-filling  analysis 
is  the  conservation  of  the  resin  mass  at  any  instant  of  time  flowing  through 
a  porous  mold  cavity.  The  fill  factor  implicitly  defines  the  amount  of  resin 
present  inside  the  mold  at  any  time  and  its  distribution.  Thus,  considering 
a  general  Eulerian  mold  domain  Q  with  mass  flux  due  to  the  incoming  resin 
at  the  gates  (as  shown  in  Figure  1),  at  any  instant  of  time,  part  of  the  mold 
cavity  is  filled  with  resin  and  the  fill  factor  (^)  in  the  filled  regions  is  greater 
than  zero  with  the  maximum  value  of  fill  factor  (T)  being  equal  to  one.  In 
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unfilled  regions,  the  resin  fill  factor  (\Er)  is  equal  to  zero.  The  mass  of  the 
resin  inside  the  mold  at  any  instant  of  time  is  then  given  by 

[  p^dQ,  (2) 

Jn 

where  p  is  the  density  of  the  polymer  resin  and  the  definition  of  'll  considers 
only  the  mold  cavity  filled  with  resin.  Since  the  mass  conservation  is  main¬ 
tained  at  any  instant  of  time,  the  mass  conservation  equation  or  continuity 
equation  for  the  resin  at  any  instant  of  time  [11,  12]  can  be  written  as  (based 
on  application  of  Reynolds  transport  equation) 


d_ 

dt 


I  p^dQ  =  —  [  ptyu  ■  n  d  dQ. 
J  n  Jd  n 


(3) 


For  a  constant  density  resin  based  on  isothermal,  incompressible  conditions, 
the  modified  mass  conservation  equation  for  the  resin  can  be  written  as 


8_ 

dt 


J 

Jn 


VdQ  = 


u  •  n  d  dQ. 


(4) 


Eq.  4  forms  the  starting  point  for  the  development  of  the  present  numerical 
methodology. 


2.2  Governing  Equations 

For  polymer  resin  flow  through  porous  media  as  in  resin  transfer  molding 
flows,  the  velocity  field  is  governed  by  Darcy’s  flow  field  approximation  and 
is  dependent  on  the  permeability  of  the  fiber  preform  and  the  resin  viscosity 
and  is  given  by 

u  =  — —  V  P,  (5) 

A* 

where  K  is  the  permeability  tensor  of  the  fiber  preform  and  is  defined  differ¬ 
ently  for  2-D  and  3-D  preform  considerations.  Using  Darcy’s  approximation 
in  the  mass  balance  (Eq.  4)  and  employing  the  Greens  theorem,  the  modi¬ 
fied  mass  balance  equation  involving  the  fill  factor  and  the  pressure  field  is 
obtained  and  given  as 
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It  should  be  noted  that  the  resin  flow  inside  a  mold  cavity  based  on 
Darcy’s  flow  approximation  is  a  pressure-driven  flow  and  the  pressure  gra¬ 
dients  are  negligible  (\/P  «  0)  in  unfilled  and  in  regions  that  are  partially 
filled  where  0  <  $  <  1.  In  filled  regions,  T  =  1,  and  the  variation  in  the 
fill  factor  are  similar  to  those  in  the  VOF  method  [7-9].  Hence,  the  modified 
mass  balance  equation  governing  the  pressure-driven  flow  field  is  given  by 


This  form  of  the  mass  balance  equation  involving  the  fill  factor  and  the  pres¬ 
sure  field  is  the  governing  equation  for  the  present  finite  element  formulations. 
A  differential  form  of  the  Eq.  7  is  given  by 


(8) 


The  boundary  conditions  for  the  mass  balance  equation  can  be  stated  as 
follows: 

dP 

At  mold  surface  -r—  =  0;  (9) 

on 

At  resin  front  P  =  0; 

At  mold  inlet  P  =  Po  (prescribed  pressure),  or 

K 

un  — - (\?P)  =  u0  (prescribed  flowrate); 

At  injection  ports  ^(t  >  0)  =  1.0. 

These  forms  of  the  mass  balance  equation  involving  the  fill  factor,  pressure 
field,  and  the  associated  boundary  conditions  are  discretized  next  using  the 
finite  element  method. 


2.3  Finite  Element  Discretization 

To  determine  the  pressure  field  and  the  resin  saturation  fill  factors,  the  mold 
cavity  is  discretized  and  modeled  using  finite  elements.  For  thin  shell  molds, 
2-D  elements  with  2-D  velocity  field  and  2-D  permeabilities  are  employed. 
For  thicker  composite  sections,  3-D  velocity  fields  as  given  by  Darcy’s  law 
and  involving  3-D  permeability  tensors  can  be  employed.  Hence,  the  present 
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formulations  are  applicable  for  both  thin  and  thick  composites.  We  now 
discretize  the  mold  geometry  by  discrete  finite  elements  and  introduce  the 
finite  element  approximations  for  both  the  pressure  and  the  fill  factor.  The 
introduction  of  finite  element  discretization  for  the  fill  factor  ^(0  <  T  <  1) 
indicates  that  the  fill  factor  has  the  same  spatial  variation  as  the  pressure  in 
each  element. 

P  =  NiPi  (10) 

*  =  NM 

In  the  aforementioned  equation,  JVj  are  the  spatial  shape  functions  at  each 
node  of  the  finite  element.  Pi  and  T;  are  the  nodal  values  of  the  pressure 
and  the  fill  factors. 

With  these  discretized  approximations  and  invoking  the  traditional 
Galerkin-weighted  residual  formulation  on  the  mass  balance  equation  given 
by  Eq.  7,  we  obtain 


d  ' 

dt  j 

j  NiNjdQ,  dq  =  [-/  vAT^v^dn  P, 

(ii) 

+  J  (vP  -njNid  dCt. 

This  discretized  system  of  equations  can  be  written  as 

CV  +  KP  =  f , 

(12) 

where 

C=  [  NTN<m, 

JQ 

(13) 

K=  [  BTKBdQ, 

Jn 

(14) 

and 

f=  [  (vP-n)NTddn, 

Jan 

(15) 

where  B  defines  the  spatial  derivatives  of  the  shape  functions  employed  for 
the  finite  element  discretization. 
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Introducing  the  finite  difference  approximation  for  the  time-derivative 
term  in  the  Eq.  12,  we  have 


"+1  _  v|/n 

At  ' 


The  discretized  finite  element  system  is  then  given  by 


C 


n+ 1 


v|/n 


+  A  tKP  =  At/. 


(16) 


(17) 


For  a  lumped  C,  the  discretized  equations  can  be  written  as 

CiW+l  -  CiW  +  A tKijPj  =  A tfi.  (18) 

Eq.  18  is  a  discrete  equation  for  both  the  pressure  and  the  fill  factor  'k 
and  is  solved  in  an  iterative  manner  as  described  next. 


2.4  Solution  Strategy 

The  fill  factors  associated  with  the  nodes  and  the  pressure  field  are  solved  in 
an  iterative  manner  based  on  the  finite  element  discrete  equations  obtained 
(given  by  Eq.  18).  The  fill  factors  at  the  beginning  of  the  simulation  are 
known,  and,  at  the  nodes  forming  the  injection  ports,  the  fill  factor  is  taken 
to  be  1.  The  iterative  procedure  for  solving  the  pressure  field  and  the  fill 
factor  until  mass  conservation  at  each  time  step  is  preserved  is  described 
next  [13]. 

1.  At  the  beginning  of  each  time  step, 


(«q)“+1  =  (*,)",  (19) 

where  m  is  the  iteration  number. 

2.  Apply  the  boundary  conditions  on  K. 

•  Boundary  conditions  based  on  injection,  vent  pressures,  and  in¬ 
jection  flow  rates  are  applied. 

•  Conditions  based  on  the  fill  factor  T,.  These  are  applied  at  nodes 
where  ^  <  1.0. 
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3.  Form 


(9i)m  =  Cu(%)n  -  Cu(%)nm+1  +  A tU  (20) 

4.  Solve  kij(Pj)m  =  ( gi)m ,  where  Kij  is  the  modified  K  after  the  appli¬ 
cation  of  boundary  conditions. 

5.  Update  the  nodal  resin  fraction  field  'Jq  using  the  modified  form  of  the 
discrete  mass  balance  equation. 

Ca^iCi  =  -  A tKijPj  +  At/,.  (21) 

-  Only  a  matrix  vector  product  and  vector  additions  are  involved 


6.  Correct  for  under  filling  or  over  filling.  Since  #  can  be  either  >  1  or  < 
0,  a  correction  is  to  be  made. 


=  maX  [( 

(22) 

7.  Continue  until  convergence, 

CW^iC1!  - 

a,  {»}£•'  <€. 

(23) 

8.  Go  to  next  time  step. 

2.5  Features  and  Advantages  of  the  Present  Method¬ 
ology 

Compared  to  the  traditional  explicit  finite  element-control  volume  method¬ 
ology  for  simulation  of  resin  transfer  mold  filling,  the  present  methodology 
has  significant  physical,  algorithmic,  and  computational  advantages.  These 
can  be  summarized  as  follows: 

•  The  method  is  based  on  a  time-dependent  mass  conservation  of  resin 
equation  that  solves  for  both  the  fill  factors  and  the  pressure  fields, 
thus  giving  a  direct  transient  problem  formulation  for  mold  filling.  The 
method  accurately  captures  the  representative  physics  while  maintain¬ 
ing  improved  numerical  and  computational  attributes. 


8 


•  The  traditional  finite  element-control  volume  methodology  requires 
identification  and  specification  of  the  control  volume  regions  associ¬ 
ated  with  nodes  of  a  given  finite-element  mesh  discretization  of  mold 
cavity.  This  may  become  cumbersome  when  higher-order  elements  are 
involved.  Further,  the  physical  volumes  associated  with  the  control  vol¬ 
ume  regions  need  to  be  computed.  The  present  methodology  only  works 
with  the  finite-element  mesh  geometry,  and  the  control  volume  regions 
associated  with  nodes  need  not  be  specified.  Complex  geometries  can 
also  be  modeled  using  higher-order  elements,  and  the  finite-element 
mesh  based  on  these  higher-order  elements  can  be  made  use  of  directly 
for  subsequent  structural  stress  simulations  in  a  concurrent  engineering 
environment. 

•  The  traditional  finite  element-control  volume  methodology  involves 
computation  of  the  velocity  fields,  normals  and  bounding  regions  of  the 
control  volume  for  the  computation  of  the  flowrates,  and  subsequent 
computations  of  fill  time  to  update  the  fill  factors  and  movements  of 
the  flow  fronts.  The  pure  finite-element  methodology  involves  only  a 
matrix-vector  product  and  vector  additions  in  the  updating  of  the  fill 
factors  and  progression  of  the  resin.  This  reduces  the  amount  of  com¬ 
putation  and  provides  a  significant  computational  advantage. 

•  The  time  increments  on  the  finite  element-control  volume  methodology 
are  restricted  so  that  no  more  than  one  complete  control  volume  region 
is  updated  at  each  time  increment.  The  present  finite-element  method¬ 
ology  does  not  have  such  a  restriction  and  computes  the  position  of  the 
flow  front  at  each  of  the  selected  discrete  time  intervals. 


3  Numerical  Examples 

A  variety  of  numerical  examples  are  presented  next  to  validate  the  present 
numerical  developments  and  demonstrate  the  applicability  to  RTM  simula¬ 
tions  involving  isothermal  conditions.  The  present  developments  are  valid 
and  can  be  extended  to  nonisothermal  conditions  where  the  temperature, 
cure,  and  non-Newtonian  effects  of  the  viscosity  are  considered.  The  nu¬ 
merical  examples  involve  2-D  triangular  and  quadrilateral  element  meshes 
involving  complex  3-D  mold  cavities.  The  methodology  is  also  directly  ap- 
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plicable  to  higher-order  elements  and  in  full  3-D  mold  geometries  involving 
3-D  elements  [13-15]. 


3.1  Circular  Plate  With  Constant  Injection  Pressure 
at  the  Center  -  Comparison  with  Analytical  Solu¬ 
tion 


A  circular  plate  with  a  hole  of  radius  i?0  at  the  center  forming  the  inlet 
gate  of  the  mold  is  considered.  The  inner  radius  of  the  plate  at  radius  Rq  is 
subjected  to  a  constant  pressure  of  P,  and  the  flow  is  driven  by  this  injection 
pressure.  The  cavity  is  of  thickness  H  and  the  permeability  of  the  fiber 
preform  is  taken  to  be  K,  with  a  porosity  of  fiber  compaction  denoted  by  <f>. 
Assuming  a  perfect  complete  radial  flow,  the  time  taken  by  the  flow  front  to 
reach  a  radius  Rf  is  given  by 


t 


rRf  rlnj^ 
KpL  WJd~z 


dr. 


(24) 


For  constant  viscosity  /xq,  the  time  to  reach  a  radius  of  Rf  is  given  by 


KP\  2  \Ra 


(25) 


A  circular  plate  model  with  the  inner  radius  Rq  =  0.0015  m,  with  perme¬ 
ability  K  =  44. 0E  —  12  m2,  porosity  $  =  0.805,  viscosity  /xo  =  0.02  poise, 
and  subjected  to  an  injection  pressure  of  69, 000  Pa  is  considered.  The  ge¬ 
ometric  mesh  employed  and  the  flow  front  contours  are  shown  in  Figure  2. 
The  flow  front,  as  computed  by  the  analytical  solution,  is  compared  with 
the  pure  finite-element  formulations  and  the  explicit  finite-element-control 
volume  formulations  as  shown  in  Figure  3  (a).  The  comparison  of  flow  front 
locations  computed  using  different  time  step  sizes  employing  the  pure  finite- 
element  formulation  are  shown  in  Figure  3  (b).  The  numerical  simulations 
are  based  on  employing  a  2-D  quadrilateral  mesh  of  600  elements,  with  a 
2-D  isotropic  permeability  tensor,  and  an  outer  radius  taken  to  be  0.1  m. 
The  results  compare  well  with  the  analytical  solution,  which  is  based  on  a 
perfect  radial  flow.  The  present  formulations  compute  the  flow  front  loca¬ 
tion  at  a  given  time  with  good  accuracy,  even  when  different  time  steps  are 
employed  to  reach  that  stage.  The  total  number  of  iterations  taken  to  com¬ 
pletely  fill  the  mold  and  the  average  number  of  iterations  per  time  step  are 
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shown  in  Figure  3  (c)  and  Figure  3  (d).  It  should  be  noted  as  the  time  step 
size  increases,  the  total  number  of  iterations  decreases,  thus  improving  the 
computational  performance. 

3.2  Circular  Plate  With  Constant  Injection  Flow  Rate 
at  the  Center  -  Comparison  with  Analytical  Pres¬ 
sure 

The  same  circular  plate  as  in  the  earlier  case  is  considered.  However,  the 
plate  is  injected  with  a  constant  volumetric  injection  Q  at  the  inlet  radius, 
forming  the  inlet  gate  of  the  mold.  Since  the  injection  is  at  a  constant 
flow  rate,  the  injection  pressure  varies  with  time  as  the  mold  gets  filled  and 
provides  a  comparison  of  the  simulated  injection  pressures  with  an  analytical 
solution.  An  expression  for  the  injection  pressure  and  the  radial  location  of 
the  flow  front  is  available  and  is  used  here.  For  this  case,  the  front  radius 
location  Rf(t)  and  the  injection  pressure  P0  at  time  t  are  given  by 


and 


Po  = 


27 xKH 


log 


ML. 

Ro 


(26) 


The  mesh  and  the  physical  properties  of  the  fiber  preform  and  the  resin 
are  the  same  as  in  the  previous  example.  The  flow  rate  Q  is  taken  to  be 
2.36e  —  04^ .  The  flow  front  contours  based  on  the  present  formulations 
with  a  time  step  size  of  1.0  sec  are  shown  in  Figure  4.  The  flow  front  loca¬ 
tion  as  computed  by  the  analytical  solution  is  compared  with  the  pure  finite 
element  formulation  and  the  explicit  finite  element-control  volume  formula¬ 
tions  as  shown  in  Figure  5  (a).  The  location  of  flow  fronts  computed  using 
different  time  step  sizes  employing  the  present  formulations  and  the  analyti¬ 
cal  solution  are  shown  in  Figure  5  (b).  The  variation  of  the  total  number  of 
iterations  and  the  average  number  of  iterations  based  on  time  steps  employed 
are  shown  in  Figures  5  (c)  and  5  (d).  Being  a  constant  flow  rate  injection, 
the  injection  pressure  changes  with  time  and  the  comparison  of  the  injection 
pressure  with  time  for  the  pure  finite  element  formulation  and  the  explicit 
finite  element-control  volume  formulation  are  shown  in  Figure  6  (a).  The 
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injection  pressure  comparisons  with  different  time  step  sizes  employed  in  the 
present  formulations  as  compared  with  the  analytical  solution  are  shown  in 
Figure  6  (b).  The  present  formulations  compute  the  given  flow  front  loca¬ 
tion  at  a  given  time  with  good  accuracy  even  when  different  time  steps  are 
employed  to  reach  that  stage. 

3.3  Flow  Fronts  Advancement  and  Time  Step  Sizes 

In  the  previous  two  examples,  the  flow  front  locations  employing  the  different 
time  step  sizes  and  the  pure  finite  element  methodology  are  shown  in  Figure 
3  (b)  for  the  case  of  constant  injection  pressure  at  the  inner  radius  and  in 
Figure  5  (b)  for  the  case  of  constant  flow  rate  injection  at  the  inlet  radius. 
The  time  data  for  these  figures  were  originally  obtained  by  determining  the 
time  step  interval  at  which  the  node  locations  along  a  certain  radial  direc¬ 
tion  got  filled.  The  resolution  of  the  time  and  the  front  location  with  this 
strategy  depends  on  the  time  step  size  and  only  indicates  that  a  given  node 
location  is  completely  filled  at  the  end  of  a  certain  discrete  time  step,  while 
the  node  physically  could  be  filled  between  two  discrete  time  steps  employed 
in  the  analysis.  Hence  the  comparison  with  the  analytical  solution  shown 
earlier  indicated  that  the  flow  front  location  accuracy  varies  with  time  step 
size.  However,  during  the  analysis,  the  computed  time  for  complete  filling 
remained  the  same  even  when  employing  larger  time  step  sizes. 

Numerical  experiments  involving  analysis  up  to  certain  intermediate  times 
also  indicated  that  the  filled  and  partially  filled  regions  up  to  that  time  were 
the  same  (within  the  bounds  of  numerical  errors)  after  the  satisfaction  of 
mass  conservation  up  to  that  time.  This  clearly  suggests  that  the  predicted 
location  of  the  flow  fronts  employing  the  present  finite  element  methodology 
are  independent  of  the  time  step  size.  It  also  suggests  that  the  location  of  the 
flow  front  was  not  correctly  resolved  with  the  strategy  employed  earlier.  In 
the  two  example  cases  involving  a  radial  flow,  the  location  of  the  front  can  be 
resolved  based  on  the  volume  filled,  as  indicated  by  the  fill  factors  obtained 
using  different  time  step  sizes  used  to  reach  that  time  level.  An  alternate 
way  of  resolving  the  flow  front  location  at  a  given  time  step  for  radial  flow  is 
described  next. 

Consider  the  computed  nodal  fill  factors  at  any  instant  of  time  t  be  Fj. 
This  time  instant  can  be  reached  employing  different  time  step  sizes  and  the 
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filled  volume  of  the  resin  V  ( t )  at  that  instant  is  given  by 


V(t)  =  FiCu  (27) 

Let  R  be  the  radial  location  of  the  resin  front  at  that  time  t.  For  the  case  of 
radial  flow,  we  then  have 


7 tH(R2  -  Rl)  =  V(t)  (28) 

from  which  the  radial  front  location  at  any  instant  of  time  can  be  computed 
based  on  mass  balance.  The  radial  front  locations  employing  the  above  strat¬ 
egy  with  different  time  step  sizes  are  as  shown  in  Figure  7  (a)  for  the  case  of 
constant  inlet  injection  pressure  and  Figure  7  (b)  for  the  case  of  constant  inlet 
injection  flow  rate.  From  the  two  figures,  it  is  clear  that  the  flow  front  loca¬ 
tion  does  not  depend  on  the  time  step  size  employed.  The  computed  resin 
front  advancement  and  impregnated  region  at  any  given  time  remains  the 
same  irrespective  of  the  time  step  size  employed  to  reach  that  stage,  though, 
higher  flow  front  progression  resolution  can  be  obtained  with  smaller  time 
step  sizes. 

3.4  Circular  Plate  With  an  Eccentric  Injection  Gate 

One  of  the  advantages  of  the  pure  finite-element  formulations  is  its  computa¬ 
tional  efficiency.  Traditional  explicit  finite  element-control  volume  method¬ 
ologies  currently  employed  in  the  resin  transfer  molding  simulations  require 
computation  of  the  velocity  field  and  the  flow  rates  and  restriction  of  time 
step  increment  sizes  based  on  the  condition  that  only  a  single  control  volume 
region  is  filled  at  each  time  step.  However,  the  present  formulations  do  not 
involve  these  computations  and,  hence,  are  computationally  efficient.  In  the 
study  of  the  earlier  two  examples  in  which  the  flow  is  radial  and  with  the  mesh 
employing  quadrilateral  elements  having  a  good  radial  symmetry,  we  found 
that  the  computational  time  for  the  traditional  explicit  finite  element-control 
volume  methodology  and  the  present  formulations  were  nearly  identical.  This 
can  be  attributed  to  the  fact  that  the  previous  two  examples  with  perfect 
radial  symmetry  allowed  for  the  simultaneous  filling  of  a  row  of  control  vol¬ 
umes  in  the  explicit  finite  element-control  volume  methodology.  Thus,  the 
nearly  equal  times  are  attributed  to  this  geometry  and  flow  symmetry.  As  a 
test,  we  consider  for  this  example  the  same  circular  plate  geometry  with  the 
injection  gate  not  at  the  center,  but  at  an  eccentricity  to  the  plate  geometry. 
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The  flow  is  not  radially  symmetric  in  this  case,  and  the  mesh  is  also  not  ra¬ 
dially  symmetric.  This  being  the  case,  the  computational  advantages  of  the 
present  developments  for  the  RTM  simulations  are  clearly  seen  from  Table 
1. 

From  Table  1  and  the  comparisons  shown  for  significantly  larger  practical 
examples  shown  later,  the  computational  advantage  of  the  present  formula¬ 
tions  is  clearly  evident.  The  total  number  of  iterations,  which  is  an  indication 
of  computational  load,  is  also  shown  and  compared  with  the  explicit  finite- 
element-control  volume  methodology  in  Table  1.  The  mesh  geometry  and 
the  flow  front  contours  based  on  a  time  step  size  At  =  1.0s  are  shown  in 
Figure  8. 

3.5  Risk  Reduction  Box 

A  complex  3-D  risk  reduction  box  geometry  manufactured  by  the  RTM  pro¬ 
cess  is  considered.  Since  the  production  and  tooling  costs  are  expensive  for 
the  RTM  process,  manufacturing  simulations  are  necessary  to  avoid  costly 
trial  runs.  With  the  presence  of  corners  and  edges,  the  fiber  preform  per¬ 
meability  varies  in  those  regions.  Simulation  tools  should  take  into  account 
these  variations  in  the  permeability  variations.  The  fill  time  and  pressure  in¬ 
formation  and  the  study  of  the  flow  front  progression  with  different  injection 
conditions  help  in  optimization  of  the  mold  and  in  process  optimization  to 
manufacture  parts  without  any  dry  spots  and  within  the  residence  life  of  the 
resin  to  allow  for  complete  filling. 

From  a  computational  perspective,  the  considered  mold  geometry  with 
geometric  and  multiple  material  regions  provides  a  base  to  study  the  compu¬ 
tational  efficiency  of  the  present  method.  A  large  finite-element  mesh  con¬ 
sisting  of  4,380  nodes  and  8,670  elements  is  considered.  The  mold  is  injected 
along  the  top  edges  of  the  risk  reduction  box  geometry.  The  mesh  geometry 
and  the  flow  fill  contours  showing  the  fill  pattern  are  shown  in  Figure  9. 
Isotropic  permeabilities  of  1.0075957e-10  m2  in  the  flat  regions  and  corner 
permeabilities,  which  are  higher  in  magnitude  by  10-100  times,  were  used  for 
the  simulations.  The  fill  contours  are  representative  of  actual  conditions  and 
an  examination  of  a  representative  box  manufactured  by  the  RTM  process. 
The  computational  time  and  the  iteration  count  of  the  pure  finite  element 
methodology  for  different  time  step  sizes  are  compared  with  the  traditional 
finite  element-control  volume  methodology  as  shown  in  Table  2. 
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3.6  Crew  Capsule  of  Composite  Armored  Vehicle 
(CAV) 

A  large  finite  element  model  representing  the  crew  capsule  of  a  composite  ar¬ 
mored  vehicle  (CAV)  is  considered  in  this  analysis.  The  finite-element  model 
employed  has  12,187  nodes  and  24,240  triangular  elements.  Two  views  of 
the  mesh  geometry  are  shown  in  Figure  10.  For  a  demonstrative  injection 
location,  the  fill  pattern  contours  are  shown  in  Figure  11.  The  explicit  finite 
element-control  volume  methodology  involved  a  total  number  of  11,007  iter¬ 
ations,  and  the  total  number  of  iterations  for  the  present  methodology  are 
significantly  smaller,  as  shown  in  Table  3.  The  computational  time  of  the 
pure  finite  element  methodology  for  different  time  step  sizes  is  compared  with 
the  explicit  finite  element-control  volume  methodology,  as  shown  in  Table  3. 

4  Concluding  Remarks 

This  report  describes  a  pure  finite  element-based  methodology  for  RTM  filling 
simulations.  The  pure  finite-element  methodology  has  significant  computa¬ 
tional  advantages  compared  to  the  traditional  explicit  finite  element-control 
volume  approaches.  For  the  mold  filling  simulation  of  the  crew  capsule  with 
12,187  nodes  and  24,240  elements,  the  pure  finite-element-based  methodology 
is  faster  by  a  factor  of  nearly  30,  even  with  a  time  step  size  of  1.0  sec.  More 
importantly,  the  proposed  methodology  inherits  good  accuracy  to  accurately 
capture  the  representative  physics  while  maintaining  improved  numerical  and 
computational  attributes.  The  present  methodology  has  no  time  step  incre¬ 
ment  restrictions  as  in  the  traditional  explicit  finite  element-control  volume 
methodology  and  works  only  with  a  finite  element  mesh  geometry,  making 
its  application  direct  and  easier.  The  mold  filling  simulation  examples  illus¬ 
trate  the  practical  applicability  of  the  method  in  manufacturing  simulation 
situations.  The  developments  are  valid  and  can  be  easily  extended  to  other 
situations  involving  pressure-driven  filling  and  porous  media  flows. 
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INJECTION  PORTS 


FLOW  FRONT 
BOUNDARY 


Figure  1:  General  partially  filled  mold  cavity 
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(a)  Mold  plate  mesh  geometry 
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(b)  Fill  contours  in  the  mold 


Figure  2:  Circular  plate  geometry:  constant  pressure  injection 
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(a)  Comparison  with  explicit:  FE/CV 


(b)  Pure  FE  comparisons 


(c)  Pure  FE:  total  number  of  itera¬ 
tions 


(d)  Pure  FE:  ave.  iterations/time 
step 


Figure  3:  Circular  plate  mold-flow  front  and  iteration  comparisons:  constant 
pressure  injection 


18 


80.00 

74.29 

68.57 
62.86 

57.14 

51.43 
45.71 
40.00 

34.29 

28.57 

22.86 

17.14 

11.43 
5.714 

.00009251 


I 


~ 

a 


■ 


Figure  4:  Flow  fill  contours:  constant  flowrate  injection 
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(a)  Comparison  with  explicit:  FE/CV 


(c)  Pure  FE:  total  number  of  itera¬ 
tions 


Figure  5:  Circular  plate  mold-flow 
flowrate  injection 
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(b)  Pure  FE  comparisons 


(d)  Pure  FE:  ave.  iterations/time 
step 


and  iteration  comparisons:  constant 


(a)  Comparison  with  explicit:  FE/CV 


(b)  Pure  FE  comparisons 


Figure  6:  Comparison  of  injection  pressures:  constant  flowrate  injection 
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Figure  7:  Comparison  of  flow  fronts 


(a)  Mold  plate  mesh  geometry 


80.00 

74.29 

68.57 
62.86 

57.14 

51.43 
45.71 
40.00 

34.29 

28.57 

22.86 

17.14 

11.43 
5.714 

.00009251 


m 


(b)  Fill  contours  in  the  mold 


Figure  8:  Circular  plate:  eccentric  injection  gates-constant  flowrate  injection 
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(a)  mesh  geometry 
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(b)  Fill  contours  in  the  mold 


Figure  9:  Risk  reduction  box  mold  filling  simulation 


24 


(a)  view  1 


(b)  view  2 


Figure  10:  Crew  capsule:  mesh  geometry 


(b)  Fill  contours  -  view  2 


Figure  11:  Crew 


capsule:  fill  contours 
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Table  1:  Comparison  of  computational  times: circular  plate  with 
an  eccentric  injection  gate 


Method 

Computational  timea 

Total  no.  of  iterations 

Explicit  FE-CV 

21.96 

551 

Pure  FE 

At  =  1.0  sec 

9.27 

266 

At  =  5.0  sec 

2.80 

76 

At  =  10.0  sec 

1.77 

46 

At  =  20.0  sec 

1.41 

32 

At  =  50.0  sec 

1.00 

23 

a  Relative  to  the  actual  computational  time  corresponding  to  At  = 
50.0  sec. 


Table  2:  Comparison  of  computational  times:risk  reduction  box 


Method 

Computational  time5 

Total  no.  of  iterations 

Explicit  FE-CV 

43.86 

3808 

Pure  FE 

At  =  5.0  sec 

7.02 

557 

At  =  10.0  sec 

4.23 

329 

At  =  20.0  sec 

2.74 

202 

At  —  40.0  sec 

1.76 

136 

At  =  85.0  sec 

1.37 

97 

At  =  170.0  sec 

1.045 

74 

At  =  340.0  sec 

1.0 

59 

b  Relative  to  the  actual  computational  time  corresponding  to  At  = 


340.0  sec. 
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Table  3:  Comparison  of  computational  times: crew  capsule  of  CAV 


Method 

Computational  timec 

Total  no.  of  iterations 

Explicit  FE-CV 

104.97 

11,007 

Pure  FE 

A t  =  1.0  sec 

3.52 

350 

At  =  3.0  sec 

1.87 

176 

At  =  7.0  sec 

1.37 

124 

At  =  16.0  sec 

1.11 

93 

At  =  32.0  sec 

1.00 

83 

c  Relative  to  the  actual  computational  time  corresponding  to  At  = 


32.0  sec. 
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